1. Read the data
data<-read.csv("/Users/xueruobing/Documents/Courses/31006 Time Series Analysis and Forecasting/Final/archive/GlobalTemperatures.csv")
2. Split train and test data in period 2000-2015
data1 <- data$LandAverageTemperature[3001:3192]
ts <- ts(data1,start=c(2000,1),frequency=12)
ts
## Jan Feb Mar Apr May Jun Jul Aug Sep Oct
## 2000 2.950 4.184 6.219 9.552 11.874 14.060 14.848 14.519 12.547 9.486
## 2001 3.336 3.720 6.208 9.245 12.271 14.110 15.161 14.427 12.736 9.935
## 2002 4.026 4.704 6.772 9.229 12.248 14.084 15.354 14.560 12.950 9.988
## 2003 3.981 4.085 6.048 9.154 12.153 14.017 14.983 14.691 12.911 10.424
## 2004 3.525 4.499 6.321 9.249 11.571 13.889 14.312 14.188 12.642 10.127
## 2005 3.808 3.920 6.544 9.618 12.226 14.476 15.190 14.510 13.217 10.601
## 2006 3.286 4.430 6.329 9.055 11.786 14.443 15.042 14.913 12.875 10.289
## 2007 4.579 4.221 6.485 9.823 12.518 14.309 15.230 14.752 12.930 10.332
## 2008 2.844 3.576 6.906 9.295 12.054 14.145 15.174 14.377 12.802 10.399
## 2009 3.687 4.094 6.086 9.367 12.112 14.201 15.231 14.655 13.153 10.136
## 2010 3.737 4.399 6.738 9.671 12.406 14.421 15.213 14.768 12.863 10.442
## 2011 3.282 3.743 6.101 9.483 11.986 14.370 15.482 15.012 12.912 10.352
## 2012 3.157 3.628 6.023 9.676 12.590 14.492 15.076 14.720 13.040 10.428
## 2013 3.685 4.222 6.261 9.044 12.195 14.568 15.003 14.742 13.154 10.256
## 2014 3.732 3.500 6.378 9.589 12.582 14.335 14.873 14.875 13.091 10.330
## 2015 3.881 4.664 6.740 9.313 12.312 14.505 15.051 14.755 12.999 10.801
## Nov Dec
## 2000 6.312 3.863
## 2001 7.319 4.507
## 2002 6.892 4.038
## 2003 6.727 5.133
## 2004 7.315 4.257
## 2005 7.423 4.878
## 2006 6.955 4.987
## 2007 7.084 4.523
## 2008 7.224 4.385
## 2009 7.031 4.310
## 2010 7.487 4.292
## 2011 6.814 4.655
## 2012 7.156 4.102
## 2013 7.424 4.724
## 2014 6.713 4.850
## 2015 7.433 5.518
library(fpp)
## Loading required package: forecast
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
## Loading required package: fma
## Loading required package: expsmooth
## Loading required package: lmtest
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
## Loading required package: tseries
library(ggplot2)
tsdisplay(ts)

train_split <- window(ts, start = c(2000,1), end = c(2012,12))
test_split <- window(ts, start = c(2013,1), end = c(2015,12))
train_split
## Jan Feb Mar Apr May Jun Jul Aug Sep Oct
## 2000 2.950 4.184 6.219 9.552 11.874 14.060 14.848 14.519 12.547 9.486
## 2001 3.336 3.720 6.208 9.245 12.271 14.110 15.161 14.427 12.736 9.935
## 2002 4.026 4.704 6.772 9.229 12.248 14.084 15.354 14.560 12.950 9.988
## 2003 3.981 4.085 6.048 9.154 12.153 14.017 14.983 14.691 12.911 10.424
## 2004 3.525 4.499 6.321 9.249 11.571 13.889 14.312 14.188 12.642 10.127
## 2005 3.808 3.920 6.544 9.618 12.226 14.476 15.190 14.510 13.217 10.601
## 2006 3.286 4.430 6.329 9.055 11.786 14.443 15.042 14.913 12.875 10.289
## 2007 4.579 4.221 6.485 9.823 12.518 14.309 15.230 14.752 12.930 10.332
## 2008 2.844 3.576 6.906 9.295 12.054 14.145 15.174 14.377 12.802 10.399
## 2009 3.687 4.094 6.086 9.367 12.112 14.201 15.231 14.655 13.153 10.136
## 2010 3.737 4.399 6.738 9.671 12.406 14.421 15.213 14.768 12.863 10.442
## 2011 3.282 3.743 6.101 9.483 11.986 14.370 15.482 15.012 12.912 10.352
## 2012 3.157 3.628 6.023 9.676 12.590 14.492 15.076 14.720 13.040 10.428
## Nov Dec
## 2000 6.312 3.863
## 2001 7.319 4.507
## 2002 6.892 4.038
## 2003 6.727 5.133
## 2004 7.315 4.257
## 2005 7.423 4.878
## 2006 6.955 4.987
## 2007 7.084 4.523
## 2008 7.224 4.385
## 2009 7.031 4.310
## 2010 7.487 4.292
## 2011 6.814 4.655
## 2012 7.156 4.102
test_split
## Jan Feb Mar Apr May Jun Jul Aug Sep Oct
## 2013 3.685 4.222 6.261 9.044 12.195 14.568 15.003 14.742 13.154 10.256
## 2014 3.732 3.500 6.378 9.589 12.582 14.335 14.873 14.875 13.091 10.330
## 2015 3.881 4.664 6.740 9.313 12.312 14.505 15.051 14.755 12.999 10.801
## Nov Dec
## 2013 7.424 4.724
## 2014 6.713 4.850
## 2015 7.433 5.518
4.ets model in 2000-2015 period to check which ets model should be used.
Modelets<-ets(train_split)
Modelets
## ETS(A,N,A)
##
## Call:
## ets(y = train_split)
##
## Smoothing parameters:
## alpha = 0.2879
## gamma = 1e-04
##
## Initial states:
## l = 9.3272
## s = -5.0565 -2.4709 0.698 3.3742 5.1066 5.5886
## 4.7139 2.6274 -0.0893 -3.1369 -5.4127 -5.9424
##
## sigma: 0.3002
##
## AIC AICc BIC
## 427.6936 431.1222 473.4415
checkresiduals(Modelets)

##
## Ljung-Box test
##
## data: Residuals from ETS(A,N,A)
## Q* = 46.711, df = 10, p-value = 1.066e-06
##
## Model df: 14. Total lags used: 24
forecastets<-forecast(Modelets,h=36)
autoplot(forecastets,ylab="Global Temperature")+autolayer(test_split,series = "Actual")

accuracy(forecastets,test_split)
## ME RMSE MAE MPE MAPE MASE
## Training set 0.003932103 0.2864278 0.2247408 -0.2483583 3.563326 0.6500357
## Test set 0.165557868 0.3378695 0.2677441 2.1889085 3.904499 0.7744174
## ACF1 Theil's U
## Training set 0.1043065 NA
## Test set 0.1683787 0.1847684
5. split train and test data from 1898.12 to 2015 period
data19 <- data$LandAverageTemperature[1788:3192]
ts19 <- ts(data19, start = c(1898,12), frequency = 12)
ts19train<-window(ts19, start = c(1898,12), end = c(2012,12))
ts19test <- window(ts19, start = c(2013,1), end = c(2015,12))
lambda19 <- BoxCox.lambda(ts19)
ts19 %>% BoxCox(lambda = lambda19)
## Jan Feb Mar Apr May Jun
## 1898
## 1899 1.7387204 1.5690767 4.6398391 10.5190694 15.4066513 18.8577897
## 1900 0.4855181 2.4971332 5.9611968 10.4196748 16.0773669 19.7000082
## 1901 1.3356978 2.4891278 6.2388126 10.9339657 15.6931586 19.9519313
## 1902 1.9543782 3.1887394 5.0342606 9.9644374 15.1671459 19.5472230
## 1903 1.7374551 3.1181152 5.1047955 9.7543308 14.9852416 18.9153774
## 1904 0.6822956 1.4659815 4.4748507 9.5046686 15.0697100 19.1614687
## 1905 1.4129414 1.0183011 4.6649720 9.6242872 15.3310368 19.3503007
## 1906 1.1865281 1.4709285 5.0087905 11.1139401 15.8080466 19.7793933
## 1907 0.8671155 1.2456895 4.5763731 9.3282847 14.2201083 18.6143824
## 1908 1.4746401 2.0641766 3.9373268 9.7543308 15.6246739 19.7000082
## 1909 0.4132364 1.9827127 4.5970198 9.2528819 14.6501161 19.0614192
## 1910 1.5105828 2.0771512 4.6221124 10.2180110 15.3734458 19.4255485
## 1911 0.7084253 1.6605135 4.2424230 10.0221939 15.0495027 19.1518441
## 1912 1.5965587 2.8423132 4.2555222 10.4299488 15.8340114 19.9752173
## 1913 1.2046011 1.5640869 4.4910090 10.0783118 14.8953653 19.1441451
## 1914 2.8722997 2.6338262 5.1228320 10.2060714 15.8915344 20.0936695
## 1915 1.4684547 2.7215187 5.4265903 11.2510225 16.2785343 19.7135573
## 1916 1.8671750 2.3100652 4.9788547 10.1191618 15.3126051 19.0229659
## 1917 1.2650837 1.5441493 4.4484315 9.7932385 14.3437554 18.9710780
## 1918 1.0265807 1.6630285 4.9265431 9.4374017 14.5204900 18.9576301
## 1919 1.3271498 2.1722252 4.8251869 10.8735229 15.3070764 19.5742819
## 1920 2.1578588 2.0460322 5.5544416 10.5396561 15.6061744 19.3368007
## 1921 2.3404274 2.5332065 5.3325647 10.5928734 15.8173188 19.9849213
## 1922 1.2857307 1.5640869 5.4463422 10.6891361 15.5433081 19.7735825
## 1923 1.7033395 1.7996071 5.0087905 9.6141678 15.0862472 19.7058148
## 1924 1.5180335 2.2193411 5.7456674 9.9882119 15.7116780 19.7948904
## 1925 1.4301822 1.9762678 5.4007802 10.7959044 15.4361790 19.3059505
## 1926 2.9776440 3.3433302 6.4888632 10.3683339 15.2370801 19.6845257
## 1927 1.5590993 2.4185862 4.7138265 10.0987329 15.5562471 19.6593717
## 1928 2.4278852 2.8137365 5.1484029 10.4556419 15.6764948 19.3522294
## 1929 1.1696901 1.3125124 5.3598325 10.1992499 15.2076266 18.9960576
## 1930 1.5953081 2.5693593 5.9472988 10.3495209 15.5580957 19.8045772
## 1931 2.0538055 1.9852915 5.5437672 10.3170409 15.4011161 20.2433864
## 1932 3.1859650 2.3523246 5.1198252 11.3623097 15.8061922 19.8743499
## 1933 1.3088562 1.9196851 5.0132835 10.2521389 15.4509469 19.3464433
## 1934 1.6894676 3.3181804 4.9369976 10.0664031 16.2505654 19.9111937
## 1935 1.3161698 3.9215283 5.8178652 9.6512823 15.1910638 19.7638986
## 1936 1.4067906 1.7273369 5.2569505 10.2555529 15.8284468 19.8355814
## 1937 1.6605135 3.1153509 4.9863357 10.2777489 16.1462307 20.3017764
## 1938 2.2482082 2.9283167 6.2076816 11.6134630 16.2095613 19.7542157
## 1939 2.2771310 2.6769416 5.0987858 10.5997439 16.3363662 20.1616957
## 1940 1.5130658 2.9365287 5.8563337 11.3396877 15.9899735 20.2628460
## 1941 2.0201525 3.2777269 5.5330964 11.0256138 16.2785343 20.3543555
## 1942 2.5038074 2.1787600 5.6782245 10.7511061 15.9899735 20.2531157
## 1943 1.5890572 3.0601478 5.4189969 11.1347422 16.2636164 19.4448520
## 1944 3.0477492 3.1679410 6.0214889 10.7821165 15.9825400 20.0587022
## 1945 1.7513809 1.9055766 5.4889279 11.2093435 15.5433081 19.7677721
## 1946 2.3391061 2.9886235 5.5681712 11.4807798 16.0215732 19.7096860
## 1947 1.8429120 2.7377571 6.4418438 11.2180242 15.9732491 20.0218051
## 1948 2.8586633 2.6607605 5.0882719 10.8942390 16.3606301 20.4206040
## 1949 2.6594128 2.0512139 5.3265086 10.6272340 15.9360960 19.5646171
## 1950 1.1145721 2.0020654 5.5849601 10.2487252 15.9138122 19.8724111
## 1951 1.2711520 1.5093414 5.3477105 10.6667748 16.1760274 19.6168178
## 1952 2.4318721 2.8996035 5.1604440 10.8286638 16.0159959 19.9888032
## 1953 2.3192995 3.2206733 6.3027186 11.5767715 16.3177064 20.3192998
## 1954 1.1757003 2.4784601 5.4341857 10.3084968 15.3420978 19.9073148
## 1955 2.9214761 2.5492646 4.5911191 10.2709184 15.7190869 20.2180942
## 1956 1.6768704 2.0797476 5.0432554 9.9067438 14.8513824 19.4429215
## 1957 1.7046014 2.5265202 5.2645033 10.6650550 15.8432866 20.6528050
## 1958 3.1443962 3.3125958 5.9010140 10.9305101 16.3270358 19.4950586
## 1959 2.3127030 2.8096579 6.5249533 11.4493973 15.9621014 20.2511698
## 1960 1.9917405 3.4568877 4.1784825 10.2657961 15.1285247 20.0703566
## 1961 2.2692375 3.2958512 6.1113547 11.1226068 16.5231910 20.4946964
## 1962 2.3708496 3.3139918 6.0199416 10.9650752 15.6006254 19.9596927
## 1963 2.2363922 3.7712684 5.3568015 10.3341333 15.4472546 19.9305906
## 1964 2.2101679 2.4039846 5.0072929 9.9067438 15.1892237 20.0218051
## 1965 2.1852978 2.4146026 5.4858841 10.0221939 15.7746753 19.7290443
## 1966 1.7071255 2.7553659 6.0153001 10.3136231 15.7339067 20.2706309
## 1967 1.9878707 2.2416426 5.8132520 10.8234900 16.5381548 19.5085801
## 1968 1.5354374 2.6944885 6.9401844 10.8183167 15.4897256 19.5278997
## 1969 1.0645275 1.6567420 5.4083691 11.0273444 16.0197141 20.0528754
## 1970 2.1513333 3.3685111 5.6231499 11.0740901 15.7450234 20.1267052
## 1971 2.3232588 2.1918385 5.2116738 10.6186419 15.7542885 19.3252307
## 1972 0.9652591 1.8582308 5.6460860 10.6152054 15.6505802 20.1033848
## 1973 2.4172582 3.5979362 6.4716156 11.6571721 16.7948157 20.7329258
## 1974 1.4215584 1.8237890 5.4980610 10.6994598 15.7302015 19.8724111
## 1975 2.4451691 2.8982374 6.2933591 11.0411911 16.7441723 20.2122585
## 1976 2.2140986 2.4079655 4.4910090 10.5413720 15.2573358 19.3001672
## 1977 1.7096501 3.0023570 6.2871209 11.5034550 16.8154568 20.8346297
## 1978 1.9814235 2.9817606 6.1392926 11.4006127 16.1313364 19.5588187
## 1979 1.9479468 2.1578588 5.9334065 10.8183167 15.7320541 20.0761843
## 1980 2.3087465 3.2526588 5.7686917 11.6519253 16.7966920 20.2608999
## 1981 3.4386204 3.7726948 7.1343660 11.6834125 16.4129145 20.6176492
## 1982 1.7932511 2.7431733 5.0192754 10.8355629 16.0085601 19.6071490
## 1983 2.9666710 3.4793925 6.2824430 11.2301793 16.7197990 20.0587022
## 1984 2.2889791 2.4904618 6.2388126 10.7786701 17.1840640 19.8821054
## 1985 2.3801203 1.9801344 5.8979306 11.1260738 16.1536789 20.2356035
## 1986 2.7892788 3.7299432 6.3698678 11.5401022 16.4428059 19.9829804
## 1987 2.5519425 3.9344535 5.5819069 11.0602356 16.3382324 20.9344722
## 1988 2.9790161 2.9351598 6.6460396 11.7516901 16.9506830 21.2718831
## 1989 1.9878707 3.2178944 6.4277500 11.2232332 15.8822537 19.9985085
## 1990 2.6917878 3.2401366 7.9255517 12.0625123 16.9112203 20.8679001
## 1991 2.6270993 3.8055262 6.0880906 12.1629355 16.5082299 21.5945806
## 1992 3.1513182 3.6446864 6.6791331 11.0204223 16.3568968 20.2492239
## 1993 2.7702795 3.3825139 6.6334405 11.2249697 16.2449728 20.1364236
## 1994 2.3986783 2.1174482 6.3323733 11.7972506 16.6354832 21.3603333
## 1995 3.0973924 4.6960510 6.6854400 11.7166665 16.4222544 21.3996676
## 1996 2.1852978 3.7598607 6.0106592 11.2058716 16.4558866 20.4537441
## 1997 2.5867945 3.4119509 6.7738499 11.5750249 16.4652313 21.3347737
## 1998 2.8395897 4.9489505 6.7912409 12.3961016 17.7796931 22.1041554
## 1999 3.1029163 4.7523801 6.0199416 11.4825237 16.6766939 21.1383618
## 2000 2.3008364 4.0063910 7.1024660 12.7403090 16.9920440 21.1952846
## 2001 2.8178161 3.3475249 7.0849321 12.1964451 17.7417584 21.2934977
## 2002 3.7798280 4.7672229 7.9939319 12.1682253 17.6981529 21.2424158
## 2003 3.7157109 3.8641699 6.8307960 12.0361116 17.5182652 21.1108926
## 2004 3.0766920 4.4645734 7.2654283 12.2035020 16.4241225 20.8600708
## 2005 3.4709504 3.6290918 7.6240404 12.8578196 17.6564629 22.0151085
## 2006 2.7499459 4.3634874 7.2782382 11.8621444 16.8267176 21.9498523
## 2007 4.5822705 4.0597663 7.5288529 13.2241238 18.2113522 21.6852640
## 2008 2.1617755 3.1471647 8.2128190 12.2847123 17.3311848 21.3622997
## 2009 3.3014312 3.8770631 6.8909997 12.4120301 17.4407401 21.4724733
## 2010 3.3713109 4.3182022 7.9385682 12.9523334 17.9981231 21.9063693
## 2011 2.7445276 3.3797126 6.9147908 12.6176783 17.2029120 21.8056329
## 2012 2.5760630 3.2192838 6.7912409 12.9612567 18.3486857 22.0467615
## 2013 3.2986410 4.0612105 7.1694861 11.8428445 17.5977504 22.1972353
## 2014 3.3643121 3.0422413 7.3567890 12.8061607 18.3334165 21.7365528
## 2015 3.5738921 4.7079000 7.9418229 12.3165183 17.8195415 22.0724862
## Jul Aug Sep Oct Nov Dec
## 1898 2.9146381
## 1899 21.3898327 20.6410851 17.4709859 12.4864143 7.7226852 2.3734978
## 1900 21.4134379 20.7172877 17.1444962 12.9434114 6.4465429 3.2568347
## 1901 21.7444455 20.6645262 17.1991421 11.7709616 6.7169907 2.8164561
## 1902 21.5315406 20.2939891 16.8511210 11.4180312 5.6262071 1.9917405
## 1903 21.0618579 19.9305906 16.2934549 11.4371975 5.9982868 2.2324556
## 1904 20.7563872 20.0218051 16.0699255 11.4738045 7.1200077 2.8695716
## 1905 21.2286671 20.2531157 16.8605088 11.8375820 7.0451120 2.8423132
## 1906 21.3170820 20.4732432 16.8004446 12.2352677 5.9751007 3.4049383
## 1907 20.6254605 19.8336433 16.5999080 12.1770427 5.5239529 2.3087465
## 1908 21.4646001 19.8607792 17.1859486 11.5296293 5.6231499 2.4611402
## 1909 21.1952846 20.5922662 17.0014470 12.1030148 7.1295793 2.6648043
## 1910 21.6261142 20.4459455 16.7272977 11.7762185 5.6935399 1.9543782
## 1911 21.2384874 20.0994986 16.5493793 12.2282073 6.6980570 2.8914081
## 1912 20.5356646 19.2365728 15.7746753 10.6667748 6.0060190 2.4584772
## 1913 21.3465697 20.5376159 16.7760554 11.9007620 6.7928223 3.8899635
## 1914 21.6635721 20.5493243 17.0691781 12.4474405 6.3761208 3.1804174
## 1915 21.9498523 20.6274135 16.9788816 11.7376783 7.0976832 3.2262323
## 1916 21.1795788 20.2570077 16.7141754 11.8586349 6.1222167 1.4894990
## 1917 21.5866987 20.0451070 16.7497979 10.9719907 6.4716156 1.4043313
## 1918 21.0814691 19.9674547 17.0052084 12.5289576 6.2575046 2.1669994
## 1919 21.6162590 20.6586654 17.2104523 12.4014107 5.9565635 2.6661525
## 1920 21.5591161 20.4868947 16.8417343 11.8025097 6.0756896 1.8838020
## 1921 21.9360149 20.2200396 16.9826421 12.2299724 6.5406560 3.2679745
## 1922 21.6241431 20.2647922 16.8060739 11.8638992 6.8276300 3.2735467
## 1923 21.0461715 20.3173526 17.0089701 12.6834086 7.3567890 3.5724785
## 1924 21.6754035 20.4011148 17.1482638 12.1223948 7.3921069 3.0216009
## 1925 21.2070655 20.6254605 17.1897179 12.1100614 7.2958586 4.3167428
## 1926 21.4075361 20.8365865 17.3708359 12.7367513 7.2318220 2.8586633
## 1927 21.7128781 20.4868947 17.7095263 13.3388844 7.3054730 2.7242240
## 1928 21.6872362 20.4654431 17.2858909 12.6230054 7.2894503 3.5625865
## 1929 20.9932468 20.4517944 16.8661420 12.5768519 7.0276041 1.8097834
## 1930 21.9933522 20.9364309 16.9920440 12.2246774 7.6725251 3.5033302
## 1931 22.2130863 21.3131509 17.7645173 13.4286731 7.1407491 3.7812549
## 1932 21.8095816 20.5688411 17.8157457 12.9719661 6.6098287 3.0794508
## 1933 21.5413882 20.6078857 16.9544423 12.1840974 6.5815148 2.5051426
## 1934 21.7523387 20.6430383 16.7160499 12.5076825 7.3872890 3.7214027
## 1935 21.7286607 20.6371787 17.2274204 12.6762995 5.8856001 2.9365287
## 1936 22.3181553 20.8033258 17.2689116 12.5715286 7.1247932 3.6305090
## 1937 21.8155049 21.0912760 18.0114356 13.0559152 7.3215022 3.0119765
## 1938 22.0269773 21.0050056 18.0380663 13.1542851 7.6014337 3.0560140
## 1939 21.8550021 20.9442658 17.1501476 12.1488314 7.0546651 5.1484029
## 1940 22.2031791 20.7075150 17.7360695 12.5715286 6.8751473 3.9344535
## 1941 22.4471462 21.1324750 16.8830436 13.0219659 6.6807097 3.3293542
## 1942 21.5532065 20.6586654 17.5050245 12.8738598 7.3840774 3.6645501
## 1943 21.9834644 20.6352256 17.4105046 13.3747860 7.1678891 4.0771029
## 1944 21.9992853 21.1108926 18.0818333 13.1972544 6.8989284 2.9515937
## 1945 21.4960964 21.6931531 17.6545684 12.7812362 6.8703929 2.5774041
## 1946 21.8826584 20.7524766 17.3500641 12.2599853 7.1279839 2.5452486
## 1947 21.7641797 20.7387906 17.5428603 13.5617721 7.4033509 3.2582268
## 1948 21.9538061 20.8737725 17.5220486 12.7972581 7.0482961 3.0298543
## 1949 21.4685366 20.7368356 17.1482638 12.7029626 6.9894328 2.9721567
## 1950 21.3524682 20.0237467 17.1388451 12.0343520 5.7763702 3.0808303
## 1951 21.7148508 21.1187403 17.8423199 12.9505489 7.1120332 4.2147920
## 1952 22.1595980 21.0069655 17.6792006 12.0413909 5.8578734 3.2624038
## 1953 21.8490766 21.1933212 17.6053241 12.8079414 6.6350151 3.7399113
## 1954 21.5512367 20.8189766 17.4142834 12.7296364 7.6806113 3.1236451
## 1955 21.5532065 21.6339990 17.1426125 13.0666396 6.7122564 2.7770626
## 1956 21.4311452 19.4989217 16.8210870 11.6396847 6.3433051 3.1776442
## 1957 21.3308419 21.2384874 17.8195415 12.6407653 7.3776549 4.2147920
## 1958 21.9399683 20.2900957 16.9375268 12.4261919 7.1678891 3.4906539
## 1959 21.9617142 20.7896332 17.1275440 12.0466706 6.3651788 3.4175627
## 1960 21.4665683 20.7348807 17.4086152 12.6425415 6.4841585 4.1147086
## 1961 21.6438558 20.7818097 17.6962575 12.3430352 7.0801515 3.0422413
## 1962 21.4803471 20.5161539 17.2632526 12.5910495 7.2046378 3.6574538
## 1963 22.0744653 21.6852640 18.2819013 13.4394555 7.6515081 3.3895189
## 1964 21.6694877 20.4478951 16.4857931 11.4720608 6.7280398 2.9845055
## 1965 21.2306310 20.6996976 16.7216736 12.6265570 6.5343741 3.4119509
## 1966 22.1992165 20.5512758 17.4804397 12.0343520 6.7059449 2.8804865
## 1967 21.7286607 20.7270612 17.3311848 13.1990453 7.0339697 3.6220070
## 1968 21.1285506 20.1480868 16.5157101 12.5644315 6.5830873 2.6877375
## 1969 21.5059408 20.6938349 17.3198592 12.6283328 7.4242409 4.0438854
## 1970 21.8550021 20.4089101 17.7531371 12.3784080 6.8292130 3.0711757
## 1971 21.3505020 20.8405002 17.4917857 12.4492115 7.2238247 3.5541112
## 1972 21.5315406 21.4134379 17.0428320 12.3748698 6.4355792 3.4287908
## 1973 21.9340383 20.8209331 17.4653140 12.8417833 7.1072493 3.5047391
## 1974 21.6438558 21.1089308 16.7479226 11.9376466 6.6208456 3.1112051
## 1975 21.7681269 20.1286488 17.5163736 12.3253560 6.5862323 3.2735467
## 1976 21.3819654 20.1228181 17.0315432 10.9149625 6.2216871 3.2930618
## 1977 22.3221221 21.2542018 16.9694810 12.0396311 7.6240404 3.3881177
## 1978 21.7937875 19.7077504 17.2858909 12.3624878 7.3118839 3.3755113
## 1979 21.3013585 20.7505214 17.6962575 12.8008190 7.2094337 4.6132535
## 1980 22.5424946 21.4488553 17.8366248 12.7367513 7.8507840 3.7655638
## 1981 22.5345459 22.0665492 17.6754106 12.4598387 7.1886559 4.4939480
## 1982 22.1516760 20.6782025 17.3028734 12.1735156 6.6933251 4.2133386
## 1983 21.8806827 21.6714596 18.5550619 12.6088009 7.9825277 3.5005126
## 1984 21.9221793 21.0461715 17.2575939 12.4315034 6.4575100 2.6136536
## 1985 21.0246064 21.1933212 17.1916027 12.3607192 6.7043672 3.4386204
## 1986 21.6537135 21.0226462 16.8849218 12.2917790 6.6271425 3.3909202
## 1987 22.8369836 20.8737725 18.0761234 12.8079414 6.8434625 4.4498986
## 1988 22.6538345 22.1675206 18.5225469 13.5041849 7.2238247 4.2657152
## 1989 22.2844433 21.3426375 17.9924183 13.1113426 7.0387445 4.3634874
## 1990 22.6259895 21.2640244 17.5977504 13.3299118 8.2357508 4.3561778
## 1991 23.0303995 21.7049877 18.5148980 13.0988228 7.4564006 3.8269637
## 1992 21.3505020 20.6196020 16.6130130 12.3872542 6.5799425 3.8498526
## 1993 22.2507413 20.7544319 17.0974148 12.8506919 6.3433051 3.7926734
## 1994 22.3756870 21.0324478 18.2075403 13.5455709 8.0200100 3.9043058
## 1995 23.2881178 22.3320396 18.3009777 13.7330107 7.9760123 3.9186570
## 1996 22.7115350 22.2665998 17.7853847 12.6922960 7.6288864 4.2104320
## 1997 22.1912917 21.5807875 18.5301965 13.7384249 8.1065057 4.3444870
## 1998 23.7369617 22.8051108 18.5072496 13.4358612 7.3985316 4.7153083
## 1999 22.7632910 21.7622061 18.6028981 13.5365719 7.7016430 4.5616354
## 2000 22.7533362 22.1001963 18.2666429 12.6230054 7.2510221 3.5484629
## 2001 23.3781519 21.9182266 18.6277825 13.4250793 8.8942182 4.4763192
## 2002 23.7650703 22.1813865 19.0383454 13.5203765 8.1898991 3.7969569
## 2003 23.0224172 22.4411896 18.9633933 14.3091862 7.9206716 5.4114053
## 2004 21.6911808 21.4468874 18.4479958 13.7709191 8.8875706 4.1118138
## 2005 23.4362101 22.0823819 19.5530207 14.6318455 9.0673804 5.0267670
## 2006 23.1402104 22.8828161 18.8942579 14.0640272 8.2931320 5.1905684
## 2007 23.5163367 22.5623688 18.9999012 14.1420263 8.5052574 4.4998269
## 2008 23.4041745 21.8194540 18.7542183 14.2637249 8.7365889 4.2977775
## 2009 23.5183406 22.3697341 19.4294089 13.7871716 8.4179848 4.1886437
## 2010 23.4822763 22.5941747 18.8712238 14.3419356 9.1742504 4.1625233
## 2011 24.0223658 23.0803008 18.9653144 14.1783332 8.0624204 4.6945702
## 2012 23.2081449 22.4987834 19.2115317 14.3164626 8.6240836 3.8885298
## 2013 23.0623339 22.5424946 19.4313391 14.0042238 9.0690484 4.7969328
## 2014 22.8031191 22.8071026 19.3098062 14.1383966 7.8979049 4.9848393
## 2015 23.1581893 22.5683318 19.1325978 14.9980895 9.0840634 6.0013794
ts19new<-diff(diff(ts19,lag = 12))
tsdisplay(ts19new)

length(ts19)
## [1] 1405
length(ts19new)
## [1] 1392
6. eacf to check which arima model should be used.
Arima(0,1,1), Arima(1,1,1), Arima(0,1,2)
library(TSA)
## Registered S3 methods overwritten by 'TSA':
## method from
## fitted.Arima forecast
## plot.Arima forecast
##
## Attaching package: 'TSA'
## The following objects are masked from 'package:stats':
##
## acf, arima
## The following object is masked from 'package:utils':
##
## tar
eacf(ts19new)
## AR/MA
## 0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 x o o o o o o o o o x x x o
## 1 x o o o o o o o o o x x x x
## 2 x x o o o o o o o o o x x x
## 3 x x o x o o o o o o o x x x
## 4 x x x x o o o o o o o x x x
## 5 x x o x x o o o o o o x x o
## 6 x x x x o o o o o o o x x x
## 7 x x x x o o o o o o o x o x
7. Arima models without moving windows
auto.arima(trainnew)
## Series: trainnew
## ARIMA(0,0,2)(1,1,1)[12] with drift
##
## Coefficients:
## ma1 ma2 sar1 sma1 drift
## 0.3482 0.2826 -0.2123 -0.8501 0.0050
## s.e. 0.0872 0.0720 0.1029 0.1285 0.0021
##
## sigma^2 estimated as 0.6038: log likelihood=-175.43
## AIC=362.86 AICc=363.48 BIC=380.68
Arima001<-Arima(train_split,order = c(0,0,1),seasonal = c(0,1,1), lambda= lambda)
Arima011<-Arima(train_split,order = c(0,1,1),seasonal = c(0,1,1), lambda= lambda)
Arima111<-Arima(train_split, order = c(1,1,1),seasonal = c(0,1,1),lambda= lambda)
Arima012<-Arima(train_split, order =c(0,1,2), seasonal = c(0,1,1),lambda= lambda)
checkresiduals(Arima001)

##
## Ljung-Box test
##
## data: Residuals from ARIMA(0,0,1)(0,1,1)[12]
## Q* = 45.034, df = 22, p-value = 0.002627
##
## Model df: 2. Total lags used: 24
checkresiduals(Arima011)

##
## Ljung-Box test
##
## data: Residuals from ARIMA(0,1,1)(0,1,1)[12]
## Q* = 45.496, df = 22, p-value = 0.002295
##
## Model df: 2. Total lags used: 24
checkresiduals(Arima111)

##
## Ljung-Box test
##
## data: Residuals from ARIMA(1,1,1)(0,1,1)[12]
## Q* = 41.492, df = 21, p-value = 0.004871
##
## Model df: 3. Total lags used: 24
checkresiduals(Arima012)

##
## Ljung-Box test
##
## data: Residuals from ARIMA(0,1,2)(0,1,1)[12]
## Q* = 46.049, df = 21, p-value = 0.001259
##
## Model df: 3. Total lags used: 24
forecastarima001<-forecast(Arima001,h=36)
forecastarima011<-forecast(Arima011,h=36)
forecastarima111<-forecast(Arima111,h=36)
forecastarima012<-forecast(Arima012,h=36)
autoplot(forecastarima001) +
autolayer(forecastarima001,series="Arima(0,0,1)(0,1,1)[12]")+
autolayer(forecastarima011,series="Arima(0,1,1)(0,1,1)[12]")+
autolayer(forecastarima111,series="Arima(1,1,1)(0,1,1)[12]") +
autolayer(forecastarima012,series="Arima(0,1,2)(0,1,1)[12]") +
autolayer(test_split,series = "Actual") +
ylab("Global Temperature")

8. Prediction of 2016
Arima111new<-Arima(ts, order = c(1,1,1),seasonal = c(0,1,1),lambda= lambda)
prediction111<-forecast(Arima111new,h=12)
autoplot(prediction111) +
ylab("Global Temperature")

prediction111$mean
## Jan Feb Mar Apr May Jun Jul
## 2016 4.169204 4.451111 6.617829 9.576192 12.338021 14.421190 15.212076
## Aug Sep Oct Nov Dec
## 2016 14.795248 13.075175 10.440646 7.281139 4.814561
accuracy(forecastarima001,test_split)
## ME RMSE MAE MPE MAPE MASE
## Training set 0.06012952 0.2977159 0.2293881 0.3349614 3.612380 0.6634774
## Test set 0.13835523 0.3325615 0.2571460 2.2033342 3.998494 0.7437636
## ACF1 Theil's U
## Training set 0.00922502 NA
## Test set 0.16521064 0.1886571
accuracy(forecastarima011,test_split)
## ME RMSE MAE MPE MAPE MASE
## Training set -0.02816987 0.3019541 0.2216357 -1.267879 3.576231 0.6410547
## Test set 0.14330501 0.3274531 0.2548695 2.103706 3.855116 0.7371794
## ACF1 Theil's U
## Training set 0.06791445 NA
## Test set 0.14836004 0.1842601
accuracy(forecastarima111,test_split)
## ME RMSE MAE MPE MAPE MASE
## Training set -0.046344957 0.2781523 0.2067364 -1.390708 3.408350 0.5979601
## Test set -0.005154685 0.2829907 0.1998137 -0.192522 2.950794 0.5779371
## ACF1 Theil's U
## Training set -0.05890966 NA
## Test set 0.12775201 0.1741388
accuracy(forecastarima012,test_split)
## ME RMSE MAE MPE MAPE MASE
## Training set -0.03730304 0.2937482 0.2196464 -1.3124479 3.591380 0.6353006
## Test set 0.03928189 0.2864386 0.2065884 0.4497322 3.043145 0.5975320
## ACF1 Theil's U
## Training set 0.06165016 NA
## Test set 0.12679527 0.1696171
autoplot(forecastarima011) +
autolayer(test_split,series = "Actual")+
ylab("Global Temperature")

autoplot(forecastarima111) +
autolayer(test_split,series = "Actual") +
ylab("Global Temperature")

autoplot(forecastarima012) +
autolayer(test_split,series = "Actual") +
ylab("Global Temperature")

9. ETS model without moving windows
fit_1 <- ets(ts19train)
fit_1
## ETS(A,N,A)
##
## Call:
## ets(y = ts19train)
##
## Smoothing parameters:
## alpha = 0.2394
## gamma = 0.0451
##
## Initial states:
## l = 8.4905
## s = -2.457 0.8372 3.4847 5.27 5.7662 4.7972
## 2.6325 -0.2666 -3.4672 -5.5231 -6.0087 -5.0652
##
## sigma: 0.3049
##
## AIC AICc BIC
## 6650.206 6650.561 6728.534
checkresiduals(fit_1)

##
## Ljung-Box test
##
## data: Residuals from ETS(A,N,A)
## Q* = 81.723, df = 10, p-value = 2.305e-13
##
## Model df: 14. Total lags used: 24
forecastets1<-forecast(fit_1,h=36)
autoplot(forecastets1,ylab="Global Temperature")+autolayer(ts19test,series = "Actual")

10. moving windows
k <- 156
n <- length(ts19new)
p <- 12
H <- 36
n-k
## [1] 1236
10.1 ETS
defaultW <- getOption("warn")
options(warn = -1)
st <- tsp(ts19)[1]+(k-2)/p # gives the start time in time units,
mae_1 <- matrix(NA,n-k-35,H)
mae_2 <- matrix(NA,n-k-35,H)
rmse_1 <- matrix(NA,n-k-35,H)
rmse_2 <- matrix(NA,n-k-35,H)
aicc_1 <- array(NA)
aicc_2 <- array(NA)
for(i in 1:(n-k-35))
{
### One Month rolling forecasting
# Expanding Window
train_1 <- window(ts19, end=st + i/p) ## Window Length: k+i
# Sliding Window - keep the training window of fixed length.
# The training set always consists of k observations.
train_2 <- window(ts19, start=st+(i-k+1)/p, end=st+i/p) ## Window Length: k
test <- window(ts19, start=st + (i+1)/p, end=st + (i+H)/p) ## Window Length: H
if (i<4) {
cat(c("*** CV", i,":","len(Expanding Window):",length(train_1), "len(Sliding Window):",length(train_2), "len(Test):",length(test),'\n' ))
cat(c("*** TRAIN - Expanding WIndow:",tsp(train_1)[1],'-',tsp(train_1)[2],'\n'))
cat(c("*** TRAIN - Sliding WIndow:",tsp(train_2)[1],'-',tsp(train_2)[2],'\n'))
cat(c("*** TEST:",tsp(test)[1],'-',tsp(test)[2],'\n'))
cat("*************************** \n \n")
}
fit_1 <- ets(train_1,model="ANA")
fcast_1 <- forecast(fit_1, h=H)
fit_2 <- ets(train_2,model = "ANA")
fcast_2 <- forecast(fit_2, h=H)
mae_1[i,1:length(test)] <- abs(fcast_1[['mean']]-test)
mae_2[i,1:length(test)] <- abs(fcast_2[['mean']]-test)
rmse_1[i,1:length(test)] <- (mean((fcast_1[['mean']]-test)^2))^0.5
rmse_2[i,1:length(test)] <- (mean((fcast_2[['mean']]-test)^2))^0.5
aicc_1[i] <- fcast_1$model$aicc
aicc_2[i] <- fcast_2$model$aicc
}
## *** CV 1 : len(Expanding Window): 156 len(Sliding Window): 156 len(Test): 36
## *** TRAIN - Expanding WIndow: 1898.91666666667 - 1911.83333333333
## *** TRAIN - Sliding WIndow: 1898.91666666667 - 1911.83333333333
## *** TEST: 1911.91666666667 - 1914.83333333333
## ***************************
##
## *** CV 2 : len(Expanding Window): 157 len(Sliding Window): 156 len(Test): 36
## *** TRAIN - Expanding WIndow: 1898.91666666667 - 1911.91666666667
## *** TRAIN - Sliding WIndow: 1899 - 1911.91666666667
## *** TEST: 1912 - 1914.91666666667
## ***************************
##
## *** CV 3 : len(Expanding Window): 158 len(Sliding Window): 156 len(Test): 36
## *** TRAIN - Expanding WIndow: 1898.91666666667 - 1912
## *** TRAIN - Sliding WIndow: 1899.08333333333 - 1912
## *** TEST: 1912.08333333333 - 1915
## ***************************
##
plot(1:36, colMeans(mae_1,na.rm=TRUE), type="l",col=1,xlab="horizon", ylab="MAE",ylim=c(0.22,0.35))
lines(1:36, colMeans(mae_2,na.rm=TRUE), type="l",col=2)
legend("topright",legend=c('ANA - Expanding Window ','ANA- Sliding Window '),col=1:4,lty=1)

mean(mae_1[1201,], na.rm = T)
## [1] 0.2172144
mean(mae_2[1201,], na.rm = T)
## [1] 0.2279714
mean(rmse_1[1201,], na.rm = T)
## [1] 0.2887426
mean(rmse_2[1201,], na.rm = T)
## [1] 0.2909526
plot(1:36, colMeans(rmse_1,na.rm=TRUE), type="l",col=1,xlab="horizon", ylab="RMSE",ylim=c(0.35,0.38))
lines(1:36, colMeans(rmse_2,na.rm=TRUE), type="l",col=2)
legend("topleft",legend=c('ANA - Expanding Window ','ANA - Sliding Window '),col=1:4,lty=1)

plot(1:1201, aicc_1, type="l",col=1,xlab="Iteration Number", ylab="AICc")
lines(1:1201, aicc_2, type="l",col=2)
legend("topleft",legend=c('ANA - Expanding Window ','ANA - Sliding Window '),col=1:4,lty=1)

10.2 ARIMA models
k <- 156
n <- length(ts19)
p <- 12
H <- 36
defaultW <- getOption("warn")
options(warn = -1)
st <- tsp(ts19)[1]+(k-2)/p # gives the start time in time units,
mae_3 <- matrix(NA,n-k,H)
mae_4 <- matrix(NA,n-k,H)
mae_5 <- matrix(NA,n-k,H)
mae_6 <- matrix(NA,n-k,H)
mae_7 <- matrix(NA,n-k,H)
mae_8 <- matrix(NA,n-k,H)
rmse_3 <- matrix(NA,n-k,H)
rmse_4 <- matrix(NA,n-k,H)
rmse_5 <- matrix(NA,n-k,H)
rmse_6 <- matrix(NA,n-k,H)
rmse_7 <- matrix(NA,n-k,H)
rmse_8 <- matrix(NA,n-k,H)
aicc_3 <- matrix(NA,n-k,H)
aicc_4 <- matrix(NA,n-k,H)
aicc_5 <- matrix(NA,n-k,H)
aicc_6 <- matrix(NA,n-k,H)
aicc_7 <- matrix(NA,n-k,H)
aicc_8 <- matrix(NA,n-k,H)
for(i in 1:(n-k))
{
### One Month rolling forecasting
# Expanding Window
train_1 <- window(ts19, end=st + i/p) ## Window Length: k+i
# Sliding Window - keep the training window of fixed length.
# The training set always consists of k observations.
train_2 <- window(ts19, start=st+(i-k+1)/p, end=st+i/p) ## Window Length: k
test <- window(ts19, start=st + (i+1)/p, end=st + (i+H)/p) ## Window Length: H
if (i<4) {
cat(c("*** CV", i,":","len(Expanding Window):",length(train_1), "len(Sliding Window):",length(train_2), "len(Test):",length(test),'\n' ))
cat(c("*** TRAIN - Expanding WIndow:",tsp(train_1)[1],'-',tsp(train_1)[2],'\n'))
cat(c("*** TRAIN - Sliding WIndow:",tsp(train_2)[1],'-',tsp(train_2)[2],'\n'))
cat(c("*** TEST:",tsp(test)[1],'-',tsp(test)[2],'\n'))
cat("*************************** \n \n")
}
fit_3 <-Arima(train_1, order=c(0,1,1), seasonal=list(order=c(0,1,1), period=p),
lambda=lambda19, method="ML")
fcast_3 <- forecast(fit_3, h=H)
fit_4 <-Arima(train_2, order=c(0,1,1), seasonal=list(order=c(0,1,1), period=p),
lambda=lambda19, method="ML")
fcast_4 <- forecast(fit_4, h=H)
fit_5 <- Arima(train_1, order=c(1,1,1), seasonal=list(order=c(0,1,1), period=p),
lambda=lambda19, method="ML")
fcast_5 <- forecast(fit_5, h=H)
fit_6 <- Arima(train_2, order=c(1,1,1), seasonal=list(order=c(0,1,1), period=p),
lambda=lambda19, method="ML")
fcast_6 <- forecast(fit_6, h=H)
fit_7 <- Arima(train_1, order=c(0,1,2), seasonal=list(order=c(0,1,1), period=p),
lambda=lambda19, method="ML")
fcast_7 <- forecast(fit_7, h=H)
fit_8 <- Arima(train_2, order=c(0,1,2), seasonal=list(order=c(0,1,1), period=p),
lambda=lambda19, method="ML")
fcast_8 <- forecast(fit_8, h=H)
mae_3[i,1:length(test)] <- abs(fcast_3[['mean']]-test)
mae_4[i,1:length(test)] <- abs(fcast_4[['mean']]-test)
mae_5[i,1:length(test)] <- abs(fcast_5[['mean']]-test)
mae_6[i,1:length(test)] <- abs(fcast_6[['mean']]-test)
mae_7[i,1:length(test)] <- abs(fcast_7[['mean']]-test)
mae_8[i,1:length(test)] <- abs(fcast_8[['mean']]-test)
rmse_3[i,1:length(test)] <- (sum((fcast_3[['mean']]-test)^2))^0.5
rmse_4[i,1:length(test)] <- (sum((fcast_4[['mean']]-test)^2))^0.5
rmse_5[i,1:length(test)] <- (sum((fcast_5[['mean']]-test)^2))^0.5
rmse_6[i,1:length(test)] <- (sum((fcast_6[['mean']]-test)^2))^0.5
rmse_7[i,1:length(test)] <- (sum((fcast_7[['mean']]-test)^2))^0.5
rmse_8[i,1:length(test)] <- (sum((fcast_8[['mean']]-test)^2))^0.5
aicc_3[i,1:length(test)] <- fcast_3$model$aicc
aicc_4[i,1:length(test)] <- fcast_4$model$aicc
aicc_5[i,1:length(test)] <- fcast_5$model$aicc
aicc_6[i,1:length(test)] <- fcast_6$model$aicc
aicc_7[i,1:length(test)] <- fcast_7$model$aicc
aicc_8[i,1:length(test)] <- fcast_8$model$aicc
}
## *** CV 1 : len(Expanding Window): 156 len(Sliding Window): 156 len(Test): 36
## *** TRAIN - Expanding WIndow: 1898.91666666667 - 1911.83333333333
## *** TRAIN - Sliding WIndow: 1898.91666666667 - 1911.83333333333
## *** TEST: 1911.91666666667 - 1914.83333333333
## ***************************
##
## *** CV 2 : len(Expanding Window): 157 len(Sliding Window): 156 len(Test): 36
## *** TRAIN - Expanding WIndow: 1898.91666666667 - 1911.91666666667
## *** TRAIN - Sliding WIndow: 1899 - 1911.91666666667
## *** TEST: 1912 - 1914.91666666667
## ***************************
##
## *** CV 3 : len(Expanding Window): 158 len(Sliding Window): 156 len(Test): 36
## *** TRAIN - Expanding WIndow: 1898.91666666667 - 1912
## *** TRAIN - Sliding WIndow: 1899.08333333333 - 1912
## *** TEST: 1912.08333333333 - 1915
## ***************************
##
Plot the MAE, RMSE and AICc for six models above.
plot(1:36, colMeans(mae_3,na.rm=TRUE), type="l",col=1,xlab="horizon", ylab="MAE",ylim=c(0.18,0.33))
lines(1:36, colMeans(mae_4,na.rm=TRUE), type="l",col=2)
lines(1:36, colMeans(mae_5,na.rm=TRUE), type="l",col=3)
lines(1:36, colMeans(mae_6,na.rm=TRUE), type="l",col=4)
lines(1:36, colMeans(mae_7,na.rm=TRUE), type="l",col=5)
lines(1:36, colMeans(mae_8,na.rm=TRUE), type="l",col=6)
legend("bottomright",legend=c('ARIMA(0,1,1)(0,1,1) - Expanding Window ','ARIMA(0,1,1)(0,1,1)- Sliding Window ','ARIMA(1,1,1)(0,1,1)- Expanding Window','ARIMA(1,1,1)(0,1,1)- Sliding Window','ARIMA(0,1,2)(0,1,1)- Expanding Window','ARIMA(0,1,2)(0,1,1)- Sliding Window'),col=1:6,lty=1)

plot(1:36, colMeans(rmse_3,na.rm=TRUE), type="l",col=1,xlab="horizon", ylab="RMSE",ylim=c(1.8,2.3))
lines(1:36, colMeans(rmse_4,na.rm=TRUE), type="l",col=2)
lines(1:36, colMeans(rmse_5,na.rm=TRUE), type="l",col=3)
lines(1:36, colMeans(rmse_6,na.rm=TRUE), type="l",col=4)
lines(1:36, colMeans(rmse_7,na.rm=TRUE), type="l",col=5)
lines(1:36, colMeans(rmse_8,na.rm=TRUE), type="l",col=6)
legend("bottomright",legend=c('ARIMA(0,1,1)(0,1,1) - Expanding Window ','ARIMA(0,1,1)(0,1,1)- Sliding Window ','ARIMA(1,1,1)(0,1,1)- Expanding Window','ARIMA(1,1,1)(0,1,1)- Sliding Window','ARIMA(0,1,2)(0,1,1)- Expanding Window','ARIMA(0,1,2)(0,1,1)- Sliding Window'),col=1:6,lty=1)

plot(1:36, colMeans(aicc_3,na.rm=TRUE), type="l",col=1,xlab="horizon", ylab="AICc",ylim=c(950,1050))
lines(1:36, colMeans(aicc_5,na.rm=TRUE), type="l",col=3)
lines(1:36, colMeans(aicc_7,na.rm=TRUE), type="l",col=5)
legend("bottomright",legend=c('ARIMA(0,1,1)(0,1,1) - Expanding Window ','ARIMA(1,1,1)(0,1,1)- Expanding Window','ARIMA(0,1,2)(0,1,1)- Expanding Window'),col=c(1,3,5),lty=1)

plot(1:36, colMeans(aicc_4,na.rm=TRUE), type="l",col=2,xlab="horizon", ylab="AICc",ylim=c(216,220))
lines(1:36, colMeans(aicc_6,na.rm=TRUE), type="l",col=4)
lines(1:36, colMeans(aicc_8,na.rm=TRUE), type="l",col=6)
legend("bottomright",legend=c('ARIMA(0,1,1)(0,1,1)- Sliding Window ','ARIMA(1,1,1)(0,1,1)- Sliding Window','ARIMA(0,1,2)(0,1,1)- Sliding Window'),col=c(2,4,6),lty=1)
